Quantitative Big Imaging

• To get an idea of what is going on
• To test a hypothesis
  • Does temperature affect bubble size?
  • Is this gene important for cell shape and thus mechanosensation in bone?
  • Does higher canal volume make bones weaker?
  • Does the granule shape affect battery life expectancy?

• What we are looking at
• What we get from the imaging modality

• What comes out of our detector / enhancement process

• Identify objects by eye
  • Count, describe qualitatively: “many little cilia on surface”, “long curly flaggelum”, “elongated nuclear structure”
• Morphometrics
  • Trace the outline of the object (or sub-structures)
  • Can calculate the area by using equal-weight-paper and employing the “cut-and-weigh” method

• For segmentation this model is:
  • there are 2 (or more) distinct components that make up the image
  • these components are distinguishable by their values (or vectors, colors, tensors, …)
  • For absorption/attenuation microscopy, Beer-Lambert Law \[ I_{detector} = I_{source}\exp(-\alpha d) \]
  • Different components have a different \( \alpha \) based on the strength of the interaction between the light and the chemical / nuclear structure of the material \[ I_{sample}(x,y) = I_{source}\exp(-\alpha(x,y) d) \] \[ \alpha = f(N,Z,\sigma,\cdots) \]

Start out with a simple image of a cross with added noise \[ I(x,y) = f(x,y) \]

The intensity can be described with a probability density function \[ P_f(x,y) \]

By examining the image and probability distribution function, we can deduce that the underyling model is a whitish phase that makes up the cross and the darkish background

Applying the threshold is a deceptively simple operation \[ I(x,y) = \begin{cases} 1, & f(x,y)\geq0.5 \\ 0, & f(x,y)<0.5 \end{cases} \]

• We can peform the same sort of analysis with this image of cells
• This time we can derive the model from the basic physics of the system
  • The field is illuminated by white light of nearly uniform brightness
  • Cells absorb light causing darker regions to appear in the image
  • Lighter regions have no cells
  • Darker regions have cells

While scalar images are easiest, it is possible for any type of image \[ I(x,y) = \vec{f}(x,y) \]

The intensity can be described with two seperate or a single joint probability density function \[ P_{\vec{f}\cdot \vec{i}}(x,y), P_{\vec{f}\cdot \vec{j}}(x,y) \]

A threshold is now more difficult to apply since there are now two distinct variables to deal with. The standard approach can be applied to both \[ I(x,y) = \begin{cases} 1, & \vec{f}_x(x,y) \geq0.25 \text{ and}\\ & \vec{f}_y(x,y) \geq0.25 \\ 0, & \text{otherwise} \end{cases} \]

This can also be shown on the joint probability distribution as

We can tune the angular acceptance by using the fact \[ \vec{x}\cdot\vec{y}=|\vec{x}| |\vec{y}| \cos(\theta_{x\rightarrow y}) \] \[ I(x,y) = \begin{cases} 1, & \cos^{-1}(\vec{f}(x,y)\cdot \vec{i}) \leq \theta^{\circ} \\ 0, & \text{otherwise} \end{cases} \]

Going beyond vector the possibilities for images are limitless. The following shows a tensor plot with the tensor represented as an ellipse. \[ I(x,y) = \hat{f}(x,y) \]

Once the variable count is above 2, individual density functions and a series of cross plots are easier to interpret than some multidimensional density hypervolume.

• Here we have a shale sample measured with X-ray tomography with three different phases inside (clay, rock, and air).
• The model is that because the chemical composition and density of each phase is different they will absorb different amounts of x-rays and appear as different brightnesses in the image

Ideally we would derive 3 values for the thresholds based on a model for the composition of each phase and how much it absorbs, but that is not always possible or practical.

• While there are 3 phases clearly visible in the image, the histogram is less telling (even after being re-scaled).

For this exercise we choose arbitrarily 3 ranges for the different phases and perform visual inspection

The relation can explicitly be written out as \[ I(x) = \begin{cases} \text{Void}, & 0 \leq x \leq 0.2 \\ \text{Clay}, & 0.2 < x \leq 0.5 \\ \text{Rock}, & 0.5 < x \end{cases} \]

Matlab

The simplist is a single threshold in Matlab:

thresh_img = gray_img > thresh

A more complicated threshold:

thresh_img = (gray_img > thresh_a) & (gray_img > thresh_b)

Python (numpy)

thresh_img = gray_img > thresh

Python

thresh_img = map(lambda gray_val: gray_val>thresh,
                gray_img)

Java

boolean[] thresh_img = new boolean[x_size*y_size*z_size];
for(int x=x_min;x<x_max;x++)
  for(int y=y_min;y<y_max;y++)
    for(int z=z_min;z<z_max;z++) {
      int offset=(z*y_size+y)*x_size+x;
      thresh_img[offset]=gray_img[offset]>thresh;
    }

In C/C++

bool* pix = malloc(x_size*y_size*z_size * sizeof (bool));

for(int x=x_min;x<x_max;x++)
  for(int y=y_min;y<y_max;y++)
    for(int z=z_min;z<z_max;z++) {
      int offset=(z*y_size+y)*x_size+x;
      thresh_img[offset]=gray_img[offset]>thresh;
    }

Partial Volume Effect

• The partial volume effect is the name for the effect of discretization on the image into pixels or voxels.
• Surfaces are complicated, voxels are simple boxes which make poor representations
• Many voxels are only partially filled, but only the voxels on the surface
• Removing the first layer alleviates issue

• Shown as a function of radius

Returning to the original image of a cross

We can now utilize information from neighborhood voxels to improve the results. These steps are called morphological operations.

Like filtering the assumption behind morphological operations are

• nearby voxels in real images are related / strongly correlated with one another
• noise and imaging artifacts are less spatially correlated.

Therefore these imaging problems can be alleviated by adjusting the balance between local and neighborhood values.

A neighborhood consists of the pixels or voxels which are of sufficient proximity to a given point. There are a number of possible definitions which largely affect the result when it is invoked.

• A large neighborhood performs operations over larger areas / volumes
  • Computationally intensive
  • Can smooth out features
• A small neighborhood performs operations over small areas / volumes
  • Computationally cheaper
  • Struggles with large noise / filling large holes

The neighborhood is important for a large number of image and other (communication, mapping, networking) processing operations:

• filtering
• morphological operations
• component labeling
• distance maps
• image correlation based tracking methods

Formally these have two effectively equivalent, but different definitions.

• Pixels which share a face (red line) with the pixel
• Pixels which share a vertex (yellow dot) with the pixel

• Pixels which are distance 1 from the pixel
• Pixels which are distance \( \sqrt{2} \) from the pixel

In 3D there is now an additional group to consider because of the extra-dimension

• Voxels which share a face (red line) with the voxel (6-adjacent)
• Voxels which share an edge (yellow dot) with the voxel (18-adjacent)
• Voxels which share a vertex (purple dot) with the voxel (26-adjacent)

• Voxels which are distance 1 from the voxel
• Voxels which are distance \( \sqrt{2} \) from the voxel
• Voxels which are distance \( \sqrt{3} \) from the voxel

Erosion

If any of the voxels in the neighborhood are 0/false than the voxel will be set to 0

• Has the effect of peeling the surface layer off of an object

Dilation

If any of the voxels in the neigbhorhood are 1/true then the voxel will be set to 1

• Has the effect of adding a layer onto an object (dunking an strawberry in chocolate, adding a coat of paint to a car)

Erosion

Taking an image of the cross at a too-high threshold, we show how dilation can be used to recover some of the missing pixels

Dilation

Taking an image of the cross at a too-low threshold, we show how erosion can be used to remove some of the extra pixels

Opening

An erosion followed by a dilation operation

• Peels a layer off and adds a layer on
• Very small objects and connections are deleted in the erosion and do not return the image is thus opened
• A cube larger than several voxels will have the exact same volume after (conservative)

Closing

A dilation followed by an erosion operation

• Adds a layer and then peels a layer off
• Objects that are very close are connected when the layer is added and they stay connected when the layer is removed thus the image is closed
• A cube larger than one voxel will have the exact same volume after (conservative)

Opening

Taking an image of the cross at a too-low threshold, we show how opening can be used to remove some of the extra pixels

Closing

Taking an image of the cross at a too-high threshold, we show how closing can be used to recover some of the missing pixels

• For many applications morphology appears to be the same as filter
• The key difference is the avalanche or non-linear nature of the results
  • A single off voxel can turn its entire neighborhood off (erosion)
  • A single on voxel can turn its entire neighborhood on (dilation)
• The effects of this are most pronounced when performed iteratively
• This is very useful for filling holes, connecting objects, creating masks

A segmented slice taken from a cortical bone sample. The dark is the calcified bone tissue and the white are holes in the image

Applying the closing with even larger windows will close everything but being to destroy the underlying structure of the mask

We can characterize this by examining the dependency of the closing size and the volume fraction (note scale)

• Be careful when using large area opening / closing operations (always check the images)
• Noise and defects in images can be amplified with these larger operations

Using a better kernel shape (circular, diamond shaped, etc) the artifacts from morphological operations can be reduced

We can characterize this by examining the dependency of the closing size and the volume fraction (note scale)

• Alternative techniques
  • Convex Hull - assuming convex shapes
  • Flood Filling - fill connected pixels like in Microsoft Paint

Courtesy of Alberto Alfonso

• Understand the progress of Alzheimer's disease as it relates to plaque formation in different regions of the brain

• First identify different regions
  • Manually contoured

• The cortex is barely visible to the human eye
• Tiny structures hint at where cortex is located

• A simple threshold is insufficient to finding the cortical structures
• Other filtering techniques are unlikely to magicially fix this problem

We see that the dilation and erosion affects are strongly related to the surface area of an object: the more surface area the larger the affect of a single dilation or erosion step.